arXiv: 1501.03892vl [astro-ph.CO] 16 Jan 2015 


PROCEEDINGS 

OF SCIENCE 



Weak Lensing Simulations for the SKA 


Prina Patej}^'^, Ian Harrison'^, Sphesihle Makhathini^^, Filipe Abdalla'"'^, David 
Bacon^, Michael L. Brown'^, Ian Heywood^^, Matt Jarvis^^, Oleg Smirnov^^ 

“University of the Western Cape, ^SKA South Africa, “Jodrell Bank, ‘^Rhodes University, 

‘^University College London, ^ICG, ^CSIRO, ^University of Oxford, 

E-mail: prina83@gmail. com 

Weak gravitational lensing is a very promising probe for cosmology. Measurements are tradition¬ 
ally made at optical wavelengths where many highly resolved galaxy images are readily available. 
However, the Square Kilometre Array (SKA) holds great promise for this type of measurement 
at radio wavelengths owing to its greatly increased sensitivity and resolution over typical radio 
surveys. The key to successful weak lensing experiments is in measuring the shapes of detected 
sources to high accuracy. In this document we describe a simulation pipeline designed to simulate 
radio images of the quality required for weak lensing, and will be typical of SKA observations. 
We provide as input, images with realistic galaxy shapes which are then simulated to produce 
images as they would have been observed with a given radio interferometer. We exploit this 
pipeline to investigate various stages of a weak lensing experiment in order to better understand 
the effects that may impact shape measurement. We first show how the proposed SKA 1-Mid ar¬ 
ray configurations perform when we compare the (known) input and output ellipticities. We then 
investigate how making small changes to these array configurations impact on this input-outut 
ellipticity comparison. We also demonstrate how alternative configurations for SKA 1-Mid that 
are smaller in extent, and with a faster survey speeds produce similar performance to those orig¬ 
inally proposed. We then show how a notional SKA configuration performs in the same shape 
measurement challenge. Finally, we describe ongoing efforts to utilise our simulation pipeline 
to address questions relating to how applicable current (mostly originating from optical data) 
shape measurement techniques are to future radio surveys. As an alternative to such image plane 
techniques, we lastly discuss a shape measurement technique based on the shapelets formalism 
that reconstructs the source shapes directly from the visibility data. We end with a discussion of 
extensions to the out current simulations and concluding remarks. 


Advancing Astrophysics with the Square Kilometre Array 
June 8-13, 2014 
Giardini Naxos, Italy 


* Speaker. 


(c) Copyright owned by the author(s) under the terms of the Creative Commons Attribution-NonCommercial-ShareAlike Licence. 


http://pos.sissa.it/ 








SKA Weak Leasing Simulations 


Prina Patel 


1. Introduction 


Weak gravitational leasing is a promising probe for cosmology. Light rays from distant sources 
are bent by the gravitational potential of objects on the path to an observer, leading to a coherent 
ellipticity or shear on images of galaxies near each other on the sky. We can measure statistics 
of this shear for galaxies in the Universe, at different redshifts and angular separations. These 
statistics are sensitive to the growth history of density fluctuations in the Universe (and therefore 
the matter power spectrum), and to the expansion history of the Universe (and hence for instance, 
dark energy parameters). 


Weak shear measurements are already maturing at optical frequencies (e.g. Kilbinger et al. 


20131, and a range of future optical experiments are planned to provide tight constraints on cos¬ 
mological parameters using this probe (for instance the ground based Large Synoptic Survey Tele¬ 


scope (LSST), ( |LSST Science Collaboration et aLj |2009[ ), and Euclid space telescope, ( [Laureijs 
et al.[ 2012| ). See also [Bacon et al. (2014l and Kitching et al. (2014i). In addition, shear can be 
measured at radio frequencies (e.g. Chang et’S^j 2004| Patel et aLj 20101, and as shown in Brown 


et al. (2014) SKA will be able to provide competitive gravitational leasing measurements. 

In order to demonstrate this, it is necessary to simulate realistic images which could be ob¬ 
tained using particular SKA measurement sets, including the effects of realistic gravitational leas¬ 
ing and telescope effects. Shear measurement techniques should be carried out on these images, 
and we should verify that we can obtain faithful estimates of the true shear. In this chapter, we 
describe our efforts to make such simulations and confirm that SKA configurations which are be¬ 
ing considered by the community are suitable for making the SKA a powerful gravitational leasing 
telescope. 


2. Simulation Pipeline 

The simulations we utilise are based on a pipeline built in two parts which we briefly describe 


here and refer to Patel et al. (20131 for further details. The pipeline works by taking input images 
that contain realistic galaxy shapes and running them through a simulator. Within the simulator we 
can define the interferometer to use and other observation details such as the frequency, bandwidth, 
integration time etc. The simulator then predicts the visibilities for the given image and observation. 
The last part of the simulation then takes these visibilities and images them as would be done with 
real observed data and a restored image is produced. Once the simulation has produced a final 
restored image we then re-analyse these images with our chosen shape measurement technique and 
compare the input and out ellipticities. 

2.1 Simulation 

We firstly describe the input images we have created that are used throughout this work. The 
input images are based on a shapelet (see 


based method that is described in detail in Re 


fregier (20031 and Refregier & Bacon (20031. Briefly, the shapelets method decomposes a galaxy 
image into a series of localised basis functions called shapelets. The shapelets are a complete and 
orthonormal set of basis functions consisting of weighted Hermite polynomials, corresponding to 
perturbations around a circular Gaussian. The process of generating these input images is described 
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in detail in ( |Rowe et’SIj 2013| ), resulting in shapelet models that represent simulated, but realistic, 
(known) galaxy shapes as would have been observed with the Hubble Space Telescope (HST). Al¬ 
though these are simulated optical galaxies we make use of them as there exists no large enough 
sample of highly resolved galaxy images in the radio. In Patel et al. (2010) it was found that on 
a case-by-case basis the intrinsic shapes of radio and optical sources were only weakly correlated, 
but that the overall distribution of ellipticities were very similar at the two wavelengths. To keep 
computation time low we have created 100 such images that are 0.85 x 0.85 arcminutes^ with a 
pixel scale of 0.05 arc seconds, that each contain ~ 100 sources each. This gives a resulting input 
number density of n ~ 140 arcminutes^^, which is far larger than any current lensing survey. 

These images are then fed into the simulator which in combination with the chosen interfer¬ 
ometer and observation particulars, predicts the visibilities. If desired we are then able to include 
affects that would effect the visibilities in a variety of ways using the Radio Interferometer Mea¬ 
surement Equation (RIME) formalism (e.g. Smirnov 20111, which relates the propagation of the 
signal from the source to detector via various observational effects. We are also able to include 
Gaussian measurement noise on the visibilities. Note in our simulations we do not currently em¬ 
ploy any observation effects but we do include Gaussian measurement noise. These visibilities are 
then imaged using standard techniques (e.g. CEEAN) that we are again able to control, generating 
the output restored image. 


2.2 Shape Measurement 


We take the restored images produced by the simulation and then analyse them using the image 
based shapelet method described in ^4.2| and compare the input and output ellipticities. Note, that 
in producing the restored image using CEEAN, a convolution is performed between the model 
image and the main lobe of the synthesised beam (PSE). We simulate this PSE as well and perform 
a deconvolution within the shapelet analysis, further details on deconvolution with shapelets can 


be found in 

Rowe et al. 

(|2013|, and further details about the image plane shapelet analysis can 

be found in 

Patel et al. 

20131. We are then left with shapelet models (i.e. shapelet coefficients 


fn,m) both before simulation and post, from which we estimate a 2-component (complex) ellipticity 
according to 


£ = 


s/2f' 


2,2 


(/o,0 - /4,o) 


( 2 . 1 ) 


( [Massey et al.[|2007] ). This ellipticity estimator is the Gaussian weighted quadrupole moment cast 
into shapelet space. We then follow [Heymans et al.|(|2006|l and fit a linear model to our data points 


Ei-el =mie\ +ci, 


( 2 . 2 ) 


where e] is the true input ellipticity and £, is the measured ellipticity. In all the follows the relative 
merit of each experiment can be compared through the calculated values of m, and c,. Eor a perfect 
experiment where we fully recover exactly all the shapes in the images m,- = 0,c,- = 0. A non-zero 
nii is indicative of a calibration bias resulting from poor correction of factors that circularise images, 
poor PSE correction for example, c; / 0 suggests a systematic where even circular objects appear 
elliptical. We have made no attempt to optimise any of our analyses to reduce the biases in anyway, 
neither have we looked for the origin of these biases. As such, the bias values presented here are 
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only meaningful in a comparative sense and should not be taken to represent the final performance 
that any such experiment might achieve. In a full analysis one would hope to understand the nature 
and origin of such biases to high precision in order to correct for them. One clear use of our 
pipeline is to asses the levels of bias that may be introduced by observation effects (e.g. Direction 
Dependant Effects (DDEs)). 

The key to making a weak lensing measurement requires accurate measurements of many 
galaxy shapes. Our simulation pipeline offers a way in which all parts of the data processing 
pipeline, from raw visibilities to restored images, can be explored. In the rest of this document we 
present investigations that have been carried out to address some of the key questions most relevant 
to weak lensing studies with the SKA. 

3. SKA Baselines Configurations 

Since the configuration for SKA 1-Mid is yet to be completely finalised we explore how chang¬ 
ing the array configuration as proposed by a small amount would effect the ability to accurately 
measure the shapes of galaxies. We generated many different SKAl-Mid array configurations and 
ran them through our pipeline. In this section we describe what impact minor changes had on the 
calculated calibration and additive biases. 

3.1 Baselines Changes and Impact on Weak Lensing 

We initially calculated the bias values for the two proposed SKAl-Mid configurations. The 
first is that proposed by the SKAO (referred to as SKAM) and the other, a short time after, by 
Robert Braun (SKAM12), both are shown in Eigure[T] Both these arrays contain 254 dishes with 
197 of them within a 4 km core, and the other 57 divided into 3 logarithmically spaced spiral arms 
extending out to 100 km. Eor all simulations run in this section we have adopted an 8 hour obser¬ 
vation at 700MHz and 10 50 MHz channels pointing at declination 5 = —40°. We add Gaussian 
measurement noise to the visibilities resulting in the sources in the output images having a signal- 
to-noise of ~ 10. We have also adopted a uniform weighting scheme through this work. Also 
shown in Eigure [T] are the recovered ellipticity distributions derived for both configurations. Un¬ 
surprisingly, both these two configurations produce similar calibration values, with m, ~ —0.261 
and very small additive bias. We use these recovered values of m,- and c, as our base values to 
which we can compare the values derived from modified SKAl-Mid configurations. Note that all 
the calculated bias values for all considered arrays are given in Table [T] 

The changes we explored were: changing the spacing in the arms, taking dishes from the core 
and redistributing them into the 3 spiral arms and adding new dishes to the spiral arms. In the 
former case, we looked at changing the arm distribution from logarithmic to equidistant and linear. 
In the case of the latter two, 9, 21, 30, 39, 51, and 60 dishes were added to the arms while keeping 
the maximum extent of the arms the same. In Table[T]we show the ellipticity recovery performance 
of these other SKAl-Mid configurations described above. The first two entries in blue and magenta 
correspond to the SKAO and Robert Braun configurations respectively. The cyan rows correspond 
to the configurations with equidistant and linear arms spacing. Entries in green and yellow are 
those where either dishes were redistributed from the core (green), or new dishes (yellow) were 
added to the spiral arms. 
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Figure 1: Left panel: Array configurations for SKAl-Mid as proposed by the SKAO (April 2013) and 
Robert Braun (September 2013). Right Panel: Recovered ellipticity (only component 1 of the ellipticity is 
shown) distributions for these SKAl-Mid arrays. 



Figure 2: Cross-sections through the PSFs obtained from configurations described in § 3. 


Since these are only modest changes to the configuration (i.e. ~ 20% movement/addition of 
dishes) we see no significant improvements in performance of the recovered ellipticity distribu¬ 
tions. Since the deconvolution of the PSF in known to be one of the major causes of systematic 
error in shear measurement, in Figure we show cross-sections of the PSFs for all the consid¬ 
ered configurations. As can be seen, the change in the PSF is small and so reaffirms the result of 
consistent calibration values in the absence of any other potential causes of noise. 

3.2 Alternative Configurations 

For weak leasing the main aspect of the baselines configuration is concerned with achieving 
high sensitivity at scales where we can measure the shapes of sources in the continuum. This 
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Table 1; SKAl-Mid ellipticity recovery performance results. The first two entries in blue and magenta 
correspond to the SKAO and Robert Braun configurations respectively. The cyan rows correspond to the 
configurations with equidistant and linear arms spacing. Entries in green and yellow are those where either 
dishes were redistributed from the core (green), or new dishes (yellow) were added to the spiral arms. 
“ quoted relative to the SKAO and Robert Braun configurations. 


Name 

^tolal 

Dishes 

Ncore 

Dishes 

^anm 

Dishes 

Arm 

Spacing 

"Sensitivity for 
SNR - 10 

Number Density 
n arcminute”^ 

lUi j Cj 

SKAM (SKAO) 

254 

197 

19 

Logarithmic 

1.0 

21.97 

-0.278±0.021 

-0.258±0.020 

0.001 ±0.005 

0.017±0.005 

SKAiVI12 (Robert Braun) 

254 

197 

19 

Logarithmic 

1.0 

30.43 

-0.227±0.015 

-0.280±0.016 

-5 X 10--^ ±0.004 
0.001 ±0.004 

SKAM12EQ 

254 

197 

19 

Equidistant 

1.25 

19.32 

-0.319±0.015 

-0.301±0.015 

-0.006 ±0.004 

0.010±0.004 

SKAM12LIN 

254 

197 

19 

Linear 

1.0 

30.69 

-0.297±0.015 

-0.282±0.016 

-0.006 ±0.004 

0.004 ±0.004 

SKAM9C 

254 

188 

22 

Logarithmic 

1.0 

32.20 

-0.292±0.019 

-0.264±0.017 

-0.006± 0.005 

0.002 ±0.005 

SKAM21C 

254 

176 

26 

Logarithmic 

1.0 

32.33 

-0.307±0.018 

-0.278±0.019 

-0.002 ±0.005 

0.001 ±0.005 

SKAM30C 

254 

167 

29 

Logarithmic 

1.0 

30.67 

-0.30S±0.015 

-0.286±0.017 

-0.006± 0.004 

0.001 ±0.005 

SKAM39C 

254 

158 

32 

Logarithmic 

1.0 

25.41 

-0.333±0.014 

-0.293±0.015 

-0.006± 0.004 

0.011 ±0.004 

SKAM51C 

254 

146 

36 

Logarithmic 

1.13 

31.21 

-0.314±0.016 

-0.285±0.015 

-0.005 ±0.005 

0.008 ±0.004 

SKAM60C 

254 

137 

39 

Logarithmic 

1.13 

25.72 

-0.31S±0.015 

-0.297±0.015 

-0.005 ±0.004 

0.005 ±0.004 

SKAM263 

263 

197 

22 

Logarithmic 

1.13 

33.52 

-0.234±0.016 

-0.264±0.016 

-0.007 ±0.004 

-0.003 ±0.004 

SKAM275 

275 

197 

26 

Logarithmic 

1.13 

29.61 

-0.256±0.017 

-0.256±0.018 

-0.004 ±0.005 

-0.009 ±0.005 

SKAM284 

284 

197 

29 

Logarithmic 

1.13 

29.79 

-0.256±0.017 

-0.268±0.017 

-0.002 ±0.005 

0.002 ±0.005 

SKAM293 

293 

197 

32 

Logarithmic 

1.25 

29.87 

-0.260±0.018 

-0.269±0.016 

-0.003 ±0.005 

-0.001 ±0.004 

SKAM305 

305 

197 

36 

Logarithmic 

1.25 

30.09 

-0.272±0.016 

-0.281±0.016 

-0.005 ±0.004 

0.002 ±0.005 

SKAMPLUS 

314 

197 

39 

Logarithmic 

1.38 

31.42 

-0.322±0.015 

-0.289±0.017 

-0.006 ±0.004 

0.005 ±0.005 


translates roughly to having significant sensitivity at scales corresponding to 0.5 - 1 arcseconds. For 
this reason increasing the the number of antennas in the spiral arms out to 70-80 km is beneficial 
for weak lensing while fhe lack of fhese baselines makes such a survey unfeasible. 

To accommodafe fhe ofher 3 cosmology science cases (cosmology wifh confinuum and HI 
galaxy surveys and infensify mapping) as well as weak lensing alfemafive configuralions have been 
proposed such fhaf fhe uv coverage is as full as possible ouf fo baselines of 70 - 80 km. In order fo 
achieve a smoofh fransifion befween fhe fhree secfions of fhe array (inner core, outer core and fhe 
spiral arms), fhe core is ‘puffed’ up slighfly while fhe fofal number of dishes in fhe (inner-i-ouler) 
core is preserved. The fwo proposed configuralions are referred fo using fhe following convention 
SKAlWf-y'AkBZ, where i refers fo fhe number of dishes moved from fhe ouler core fo fhe spiral 
arms, 7 is fhe number of new dishes added fo fhe arms, k is fhe maximum exlenl of fhe spiral arms 
and I fhe maximum baseline. The latter fwo are fixed in bolh cases fo k= 72km and 1= 120km. We 
shall refer fo fhe fwo configurations, SKA1W9-0A72B120 and SKA1W9-12A72B120, as A and B 
respectively. 
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Figure 3: SKAl-Mid configurations proposed by Robert Braun (blue) and also one that is smaller in extent 
and takes into consideration the site geography (red). Results for both these conhgurations are compared in 
Table 2. 


To assess the capabilities of these alternative configurations we run a new set of simulations 
at 600 MHz, 800 MHz and 1000 MHz with 1 50 MHz channel observing at declination 5 = —30°. 
To compare we also run the same simulation using the Robert Braun SKAl-Mid configuration 
discussed above, the results are presented in Table 

We notice immediately that as we increase the frequency of the observation, the calibration 
values decrease. This is because the resolution is oc A/D, where A is the wavelength of observation 
and D the maximum baseline. In this case, our PSF is effectively becoming smaller, and so the 


higher the frequency, the more accurate our shape measurement. Also, as described in Makhathini 


|et al. (in preparation! these arrays have optimised the distribution of dishes in the spiral arms such 


that the sensitivity at angular scales of 0.4 — 1 arc second at 650 MHz can be enhances without 
significantly compromising the larger scales, so we expect them to perform best at the lower part 
of the frequency space explored here. 

Also shown in Table is the same simulation ran for the so-called SKA1V8 configuration. 
This is a slightly altered configuration of SKAl-Mid that also has a smaller extent than the two 
originally proposed configurations (maximum baseline of ~ 1.50km opposed to ~ 170km), it also 
takes into account the geography of the site. This configuration is plotted (along with SKAM12 
from above) in Figure Encouragingly, even with site topology incorporated this configuration 
produces similar calibration values to the one originally proposed that did not take this into account, 
while also bringing down the maximum extent of the spiral arms. 


3.3 SKA Capabilities 

In this Section we compute the performance of SKA and compare it to SKAl-Mid. In order to 
do this we have created a SKA configuration that consists of 5 spiral arms extending out too 150 km 
(each spiral arm has 50 dishes logarithmically spaced), but we have neglected all the dishes within 
the 1 km core. This is due to the very large number of baselines involved in a full SKA simulation 
(and hence massive computation time). Instead, since weak leasing is primarily concerned with 
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Table 2; Weak lensing simulations results for SKAM12 proposed by Robert Braun and also 2 alternatives 
(SKA1W9-0A72B120 and SKA1W9-12A72B120) that give a fuller uv coverage to baselines between 70 - 
80 km. SKA1V8 is the SKAl-Mid configuration with a smaller extent and with site geography considered. 


Array Configuration 

600 MHz 

800 MHz 

1000 MHz 


m, 

Ci 

nii 

Ci 

mi 

Ci 

SKAl REF2 (Robert Braun) 

-0.560 ±0.039 

0.028 ±0.011 

-0.508±0.031 

0.052 ±0.008 

-0.424 ± 0.021 

0.021 ±0.005 

-0.491 ±0.040 

-0.004 ±0.010 

-0.434 ±0.033 

-0.007 ±0.008 

-0.400 ± 0.021 

-0.001 ±0.005 

SKA1W9-0A72B120 

-0.655 ±0.027 

0.033 ±0.007 

-0.604 ±0.020 

0.034 ±0.005 

-0.533 ±0.020 

0.016± 0.005 

-0.639 ±0.025 

-0.006 ±0.006 

-0.557 ±0.021 

-0.006 ±0.005 

-0.530± 0.021 

0.003 ± 0.005 

SKA1W9-12A72B120 

-0.582 ±0.043 

0.038 ±0.011 

-0.563 ±0.020 

0.008 ±0.005 

-0.530±0.018 

0.002 ±0.004 

-0.596 ±0.040 

-0.006 ±0.010 

-0.532 ±0.020 

-0.007 ±0.005 

-0.530± 0.017 

-0.008 ±0.004 

SKA1V8 

-0.545 ±0.042 

0.026 ±0.011 

-0.519±0.037 

0.071 ±0.010 

-0.458 ±0.023 

0.052 ±0.006 

-0.506 ±0.046 

0.016±0.011 

-0.480 ±0.036 

-0.003 ± 0.009 

-0.428 ±0.024 

0.004 ±0.020 


the longer baselines anyway, and also at this stage we are only interesting gaining some idea as to 
what SKA might be able to achieve, we adopt this simplification. The resultant SKA configuration 
is shown in the left hand panel of Figure]^ 

To compare appropriately to the SKAl-Mid simulations discussed in ^we again run the SKA 
simulation with an 8 hour observation at 700MHz and 10 50 MHz channels pointing at declination 
5 = —40°. We have also again corrupted the visibilities by an amount that results in sources being 
at a SNR ~ 10. The received ellipticity distribution is shown in the right hand panel of Figure]^ 
the resulting calibration values are: 

mi = -0.357 ±0.005 
Cl = 3x lO^'^TO.OOl 
m2 = -0.354 ±0.005 

C2 = -0.002 ±0.001. 

(3.1) 

We see that the multiplicative bias values recovered from SKA seem to be worse than for 
SKAl-Mid. We note at this stage that this is only a notional SKA configuration that we have 
simulated and so can not be completely relied upon when comparing to the more sophisticated 
SKAl-Mid configurations. Since our SKA configuration has some dishes missing, we are invari¬ 
ably missing many short and intermediate length baselines, that also carry shape information of 
scales that are relevant. We also see that error bars on m, are a factor of 10 smaller than for SKAl- 
Mid. This is due to the many more sources that reach our final catalogue as SNR ~ 10 sources, 
meaning this is a more precise measurement. 


3.4 Calibration Requirements for SKA 

To provide some context for the obtained calibration values, we calculate the requirements 
on m and c for stages of SKAl-Mid and SKA and also for comparison, current and future optical 
surveys such as the Dark Energy Survey ( |The Dark Energy Survey Collaboration 20051, Euclid and 
ESST. We adopt the requirements as computed in Amara & Refregier ( |2008| ), which are based upon 
the parameters: sky area, galaxy median redshift and galaxy number density. The requirements are 
set such that the statistical error is equal to the systematic error and thus provides an upper limit on 
the level of bias allowed. 
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Figure 4: Left panel: A mock SKA configuration with 5 spiral arms extending to 150 km, but no dishes 
within the inner 1km of the core. Right panel: Recovered ellipticity distribution for SKA. 


In addition, we follow the convention of converting the multiplicative and additive biases into 
a single quality factor Q, computed here as in Voigt & Bridle] ( 20T0) | with an assumed rms cosmic 
shear of Oy = 0.03. In Table we show the requirements for notional surveys conducted with 
an early phase of SKAl-Mid (SKAl-Mid early), SKAl-Mid and SKA, along with corresponding 
numbers for DES-like and Euclid/LSST-like surveys for comparison. 

The values quoted for the number density and median redshifts are derived for an envisaged 2 
year (net) continuum survey over 3 possible survey areas. The specification used are those given 
in Braun (2013 1 , which in turn use the SKAl-Mid baseline design and the SKADS simulations 


of Wilman et al. (2008 1 . The sensitivity levels have been chosen appropriately for weak lensing 
angular scales of 0.5 arc seconds at Band 2 and the galaxy number densities correspond to > 10a 
detections. SKAl-Mid early is defined to be such that it has 50% of the sensitivity of SKAl-Mid. 

We note how these requirements are orders of magnitude smaller than those derived in the 
preceding section. In our simulations we have not attempted to optimise any of the parameters 
(either in the simulation or the shape measurement analysis) to seek out the smallest calibration 
values, e.g. we have made no attempt to optimise the imaging of the simulated data by investigating 
other imaging methods other than CEEAN. The requirements quoted represent the levels that need 
to be achieved in order for the error budget to be equal between the systematics and statistics. 
We hope that we can utilise our pipeline further to understand the various systematics and explore 
different imaging techniques etc. to provide more robust values of for the calibration biases. 


4. Shear Measurement Techniques 

As discussed above, the signal in weak gravitational lensing is the small shearing of galaxy 
images by foreground matter. The smallness of this shearing (typically of order 1%) and its sensi¬ 
tivity to change in cosmological parameters (typically of order 0.01% for a 1% change in the dark 
energy equation of state w) means that any effect which is capable of biasing results must be care¬ 
fully controlled. One place in which such a bias may enter is in the translation from real, noisy data 
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Table 3; Requirements on multiplicative and additive biases on ellipticity measurement for proposed SKA 
weak lensing surveys to be dominated by statistical rather than systematics uncertainties, and for DES-like 
and Euclid/LSST-like for comparison. Q is calculated from m and c as in|Voigt & Bridl^(|2010|l. 


Experimenl 

^sky 

^gal 

Zm 

m < 

c < 

Q> 

DES-like 

5000 

12 

0.8 

0.004 

0.0006 

260 

Euclid/ESST-like 

20000 

35 

0.9 

0.001 

0.0005 

990 

SKAl-Mid early 

1000 

3.0 

1.0 

0.014 

0.0012 

62 

SKAl-Mid early 

5000 

1.2 

0.8 

0.012 

0.0011 

79 

SKAl-Mid early 

30940 

0.35 

0.5 

0.011 

0.0011 

80 

SKAl-Mid 

1000 

6.1 

1.2 

0.0090 

0.00095 

103 

SKAl-Mid 

5000 

2.7 

1.0 

0.0067 

0.00082 

140 

SKAl-Mid 

30940 

0.9 

0.7 

0.0058 

0.00076 

164 

SKA 

1000 

37 

1.6 

0.0031 

0.00055 

318 

SKA 

5000 

23 

1.4 

0.0019 

0.00043 

523 

SKA 

30940 

10 

1.3 

0.0012 

0.00035 

825 


to a map of shear measurements across the field-of-view. Typically this is done by measuring the 
ellipticity of galaxies identified in the data, which is changed by shear. The preponderance of opti¬ 
cal data of the quality necessary for weak lensing has led to the development of a large number of 
different techniques for performing this shape measurement process which take imageplane data as 
their inputs. Among the first derived (and subsequently most widely used) are methods which use 
weighted quadropole moment of combined galaxy-PSF images to measure ellipticities directly in a 
non-parametric way (KSB Kaiser et’ahl l 19951; KSB-i- Hoekstra et al. ( 1998| l; Re-Gaussianization 
Hirata & Seljak|((2^003|l). Another popular approach is to assert that the galaxy images may be mod¬ 


elled with some analytic brightness distribution (such as a Gaussian or Sersic profile) and find fhe 


besf fiffing paramefers, including elliplicify paramefers, for each source (IM3SHAPE Zunfz el al. 
( |2013| ), lensfif [Miller ef al.| ( |2007| )). 

In fhe radio regime fhe approach which has found mosl application is fhal of shapelels, which 
reconslrucls fhe dala using an orfhonormal sel of basis functions. How fhese basis funcfions Irans- 
form wifh shear is known, meaning fhe besl-fifling coefficienls for an image may be used lo form 
an unbiased estimator for fhe shearing if has undergone fhrough comparison wifh coefficenls from 
some ‘unlensed’ sample. Shapelels also have fhe advanfage of having similarly simple and analyfic 
Fourier fransformafions which also remain localised, facilifaling fheir use in direclly modelling 
visibilify dafa ralher lhan reconslrucfed images. 

As demonslrafed above, using simulalions wifh known elliplicify dislribulions provides a way 
of probing differenl aspecls of fhe weak lensing pipeline. Mosl nolably in fhe oplical communily 
simulalions have been used for lesfing shape measuremenl lechniques. Over fhe previous 10 years, 
fhe Shear TEsling Programme (STEP) and GRavifafional lEnsing Accuracy Tesling (GREAT) (see 
Mandelbaum el al. (2014|) and references Iherein for a brief history) iniliafives have simulaled 


large oplical weak lensing dalasels and invifed parlicipanls fo (blindly) measure fhe shear in Ihose 
images. This has allowed relafive proficiency of differenl shear measuremenl melhods and how 
Ihey read lo changes in dala parameters, such as source size, S/N and simulated galaxy model 
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complexity, to be quantified. These investigations have given insight into the behaviour of shear 
measurement algorithms in optical data, but we may expect the challenges of radio data to be 
significantly different. 

The production of images from interferometer data via deconvolution using algorithms such 
as CLEAN is a highly non-linear process and has the potential to produce spurious cosmic shear 
signals in addition to (and in convolution with) any introduced by a particular shape measurement 
algorithm. The noise in radio images is highly correlated (though optical methods have experience 
with this in dealing with multi-epoch data). The ability of presently available techniques to deal 
with these challenges is currently unclear, motivating a systematic evaluation along the model of 
STEP and GREAT programmes. A direct follow-on from these efforts may evaluate the ability of 
optical shape measurement algorithms to be extended to radio images, but a worthwhile comparison 
may also be done of the image reconstruction algorithms themselves, either in separation or in 
conjunction with shear estimation; in addition to long-established algorithms such as CEEAN and 


Maximum Entropy approaches, several new methods are under active research (e.g. Sutter et al. 
( |2014| l, |Carrillo et aLl ( |2014| l and references therein). 

If image-reconstruction algorithms are unable to perform to the low level of systematics nec¬ 
essary for weak lensing science, it may be necessary to perform shear measurement directly using 
the visibilities produced by the telescope. Indeed, the only current detection of a shear signal in 
radio sources was performed in the visibility plane with shapelets. Notionally, being able to avoid 
the deconvolution process necessary for imaging is an important strength of visibility-plane shear 
measurement methods. The deterministic nature of instrumental effects in radio astronomy is often 
touted as a key advantage of radio weak lensing and it manifests here through the ability to forward 
convolve sky models in a well-defined way when searching for the best fitting parameters. How¬ 
ever, the computational challenges of such a procedure are potentially great. Information from each 
individual source is no longer localised in the visibility plane, meaning all sources must be fitted 
simultaneously, rather than simply taking cut-outs of images around a single source as is done in 
image-plane methods. Starting with a naive model-fitting approach, if we attempt to fit 5 param¬ 
eters per source (as is typically done in mainstream optical methods) for large numbers of sources 
over (extremely) large numbers of visibilities, this quickly becomes computationally unfeasible. 
There are various potential ways of alleviating these problems: data volumes could be reduced 
with averaging of visibilities on a grid, the number of necessary simultaneous fits could be reduced 
by averaging around a single source at the phase centre or by employing methods such as visibility¬ 
stacking ( Eindroos et al.[ 20141 and we may also expect other, novel techniques to be developed. 
In addition to the computational problems involved in the analysis, the sheer volume of storage 
necessary to maintain access to unaveraged visibilities from the SKA may also be prohibitive. 


4.1 radioGREAT 

In order to investigate the issues discussed above, the radioGREAlj^ programme has been 
initiated, with three key goals: 

• What are the requirements on shape measurement for cosmology with weak gravitational 
lensing in the radio band? 

* http ://raciiogreat. j b. man. ac .uk 
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• Can we make images of the necessary fidelity to measure shapes of radio star-forming galax¬ 
ies to the level of these requirements? 

• Can we measure shapes of radio star forming galaxies to the level of these requirements 
whilst leaving our data in the visibility plane? 

For an initial challenge, the task should be kept as simple as possible, with complications intro¬ 
duced individually in order to evaluate their effect on shape measurement. We may expect to 
follow closely the structure of the GREAT08 optical/NIR challenge, with a fiducial branch con¬ 
sisting of non-overlapping, simple Sersic galaxy models with identical, high signal-to-noise values 
(> 100) and constant radii. Other branches may then consist of datasets containing individual 

modifications, to systematically evaluate the effect of e.g. altering galaxy size, reducing signal-to- 
noise, altering dynamic range and altering bandwidth and time smearing. 


4.2 Shapelets in Real and Fourier Space 

In contrast to traditional telescopes interferometers do not provide a direct image of the observed 
sky, but instead measure its Fourier transform at a finite number of uv points that correspond to 
each antenna pair in the array. The real space image must then be reconstructed from this dis¬ 
cretely sampled visibility data, while also deconvolving the effective beam that arises from the 
finite sampling. Several methods exist to perform this task, e.g. CLEAN, MEM. These methods 
are well tested and are appropriate for various applications; however, these methods are non-linear 
and do not necessarily converge in a well defined manner. 

The Hermife polynomials fhaf form fhe shapelef basis sef have some remarkable properties 
which greafly facilifale fhe modelling of fhe source shapes. Of particular inferesf here is fhe prop- 
erfy fhaf fhey are invarianf (up fo a rescaling) under of Eourier fransform and fhus are nafurally 
suifed for inferferomefric imaging. In opfical surveys, shapelefs have been liffle ufilised after some 
concepfual concerns were raised wifh fhe approach (of shapelefs like mefhods), however, similar 
concepfual problems were also raised wifh fhe KSB mefhod alfhough fhis mefhod confinued fo be 
popular in fhe liferafure, see Melchior el al. (20111 for more defails. 

In fhe real space (or image plane) applicafion of fhe shapelefs technique fhe surface brighlness 
/(x) of an objecl is decomposed as 


/(x) =£/nRn(x;i3); 

n 


(4.1) 


where ^ 

[2("'+"2)7rj32ni !n2!] ^ 


(4.2) 


are fhe two-dimensional orthonormal basis functions with a characteristic scale j3, Hm{f]) is the 
Hermite polynomial of order m, x = {xi,X 2 ) and n = (ni,n2)- The basis functions are complete 
and if j3 and x are chosen to be close to the size and location of the galaxy, then the expansion will 
yield a quick convergence. 

The Eourier transform of an objects intensity can be written as 


1 r°° 

m = { 2 n )-2 / 

J —oo 


(4.3) 
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and decomposed as 

/(k)=£/„B„(k;^), 

n 

where Bn(k; j3) obey the dual property 

B„(k;^) = /(”'+”^)B„(k;^). 


(4.4) 


(4.5) 


The invariance (unto a rescaling) under Fourier transform makes this basis set a natural choice for 
interferometric imaging. 

As mentioned above, an interferometer correlates the signals measured by antenna pairs into a 
complex visibility. Each of these visibilities occupies a point on the uv plane which corresponds to 
the projected baseline formed between the antenna pair. In practice, the visibilities are not exactly 
a two-dimensional Fourier transform of the sky brightness. The visibility measured for an antenna 
pair {i,j) at a time and frequency of t and v respectively is given by 



where /(I, v,t) is the surface brightness of the sky at location 1 = {£,m) with respect to the phase 
centre and A{i,v) is the frequency dependent primary beam. The (n, v, w) coordinates are given by 


^ M^ 


' sin//o cos//o 0 \ 

V 

= 

— sin^cos//o sin^sin//o cosSq 



k cosSqCosHq —cosSosinHo sin^ / 


1 

X 



(4.7) 


where X is the wavelength of the observation, Hq is the hour angle and 5o the declination. y ^ 
the coordinate differences between the two antennas measured in a fixed-Earth coordinate system 
in which the sky rotates about the axis. We also note that locus that is traced by the {u,v,w) 
coordinates also define the synthesised beam pattern (or PSF). Since these coordinates are entirely 
determined from the antenna and source positions and the time and frequency of the observation, 
the synthesised beam is entirely known for interferometers. Only for observations at the zenith 
(z ~ 0), in the absence of a primary beam (A ~ 1) and for small displacements from the phase 
center {£,m « 1) does the relationship between the measured visibilities and the desired sky 
brightness reduce to an extact Fourier transform. 

Since visibility datasets can be large (> 10^), directly fitting the shape parameters to all the 
uv data can be very computationally expensive. Instead one can apply a binning scheme to reduce 
the size of the data but without losing any information. In the uv plane we can grid the data using 
a cell size Am = 0.5A£^^, and average the data in each cell and similarly for Av. Ai is chosen to 
be one half of the intended field-of-view and the 0.5 factor accounts for the Nyquist frequency. 
This choice of Am, v is designed to minimise the number of cells but also to avoid smearing at large 
scales which mimic primary beam attenuation. 

In our implementation of this technique we model the intensity fs{i,m) of each source s as 
a sum of the shapelet basis functions Bn(l —centred on the source centroid ls{£s,>ns) with 
scale fig, by estimating the shapelet coefficients f„s of the sources given the binned visibility data 
described above. In principle, a fully sampled uv plane provides all the shape information about 
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Figure 5: uv coverage of the eMERLIN experiment used to demonstrate shape measurement using shapelet 
fitting in the uv plane. 



the sources and a linear decomposition similar to the image domain described above could be 
performed. However, given the finite sampling in the uv plane this is not possible here. To alleviate 
this problem we make a linear fit to the uv plane with the shapelet coefficients as free parameters. 

4.3 uv Plane shapelet Fitting 

To illustrate the technique we have performed simple simulations. We have taken a sample 
of the source models used in the studies discussed previously and made very simple sky patches 
containing a grid of 10 x 10 sources. In order to keep the simulations small (i.e. the number of visi¬ 
bilities) we have adopted the eMERLIN array for this experiment. The observation configuration is 
as follows: the fields are all observed at 5 = +60°, at 1.4 GHz with a single 125 kHz channel, and 
for 24 hours and 20 second integration. The resulting uv coverage on which the following results 
are based is shown in Figure]^ We produce ~ 50 tiles that have a grid of 100 sources each, each 
tile consists of sources with the same Umax, the range of Umax is 0 — 20. We simulate the visibilities 
for each of these sources, bin and fit for the shapelet coefficients as described above. As in the 
previous Sections we then compare the input and output ellipticities. Our current implementation 
of this technique fits all source with the same number of coefficients (i.e. for the same Umax), hence 
our motivation to tile sources according to Umax- In practice, we would want to fit all sources based 
on some guess or estimate for the what the Umax for each source in our field should be. Likewise, we 
currently input the source positions and jS parameters into the code when we do the shape fitting. 
We are working towards an implementation where the Umax and jS parameters are also fir for. 

As an initial study, in Figure we show the comparison of input and output ellipticities having 
fit all tiles with an Umax = 4, this is the lowest Umax we can chose for our adopted ellipticity estimator. 
We see that the bias we achieve in this (simplified) scenario is much less than anything we achieved 
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Figure 6: Ellipticity comparison based on shapelet fitting in the uv plane. The solid curves are show the 
binned data and the dashed curves are the best linear fits. The red lines are for ei and the blue for £ 2 . 


in the image plane analyses of the previous Sections, the calculated bias values are 

mi = 0.140 ±0.012 (4.8) 

Cl = -0.003 ±0.002 
m 2 = 0.176±0.010 
C2 = -0.014 ±0.002 


Once again we find very little additive bias while a much smaller multiplicative bias than we have 
encountered elsewhere in this chapter. Although this in not a meaningful comparison as the input 
images and simulations are fundamentally different in both cases. We aim to expand these uv plane 
simulations considerably in order to make much better comparison between the performance of 
shapelet fitting in the uv and image planes. Ultimately we hope to extend this to real data that has 
already been analysed in the the image plane in|Patel et al.|(|2010|l. 


5. Discussion/Conclusion 

Over recent times weak leasing has emerged as a very promising tool for cosmology. Making 
weak leasing measurements requires the measurement of very many galaxy shapes to high pre¬ 
cision, while also carefully controlling systematics. Due to the involved nature of making radio 
images making use of simulated data provides a important way to probe such systematic effects in 
the system while also being key to testing current and new techniques. 

In this chapter we have provided details of the ongoing efforts carried out to address the fea¬ 
sibility of weak leasing experiments with SKA 1-Mid and SKA in the context of practical shape 
measurement. We have implemented a pipeline that allows us to explore many aspects of the data 
gathering and reduction processes, allowing us to asses if these are able to meet the high demands 
on data quality that weak leasing measurements require. 

We have explored how small changes to the proposed SKA 1-Mid array configuration has little 
impact on ones ability to accurately measure the shapes of galaxies. We have introduced alternative 
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configurations that are smaller in overall extent that give similar performance to those originally 
proposed and which have a better survey speed. Also, demonstrated was the performance of SKAl- 
Mid with the site geography taken into consideration and this provides encouraging results. 

Using a notional SKA configuration we find slightly diminished calibration biases in compar¬ 
ison to SKAl-Mid but this result requires further investigation as it is unlikely this can truly be the 
case. A more complete configuration without missing baselines should give a much clearer idea of 
how SKA will perform in comparison on SKAl-Mid. 

In our present studies we have not tried to optimise the simulations or shape measurement 
analysis in order to achieve the optimal calibration biases. However, we did present the require¬ 
ments that will be required for constraining cosmology with real SKA surveys. We can make use 
of the pipeline that we have developed to determine where these biases originate. 

We introduced the forthcoming radioGREAT challenge that will investigate current and new 
methods of shape measurement and their applicability to radio data. Starting with simple simula¬ 
tions we hope to emulate the successes of the STEP and GREAT challenges in the optical to address 
the difficulties that radio data gathering and reduction will most likely pose. We also demonstrate 
the shapelet shape measurement technique as applied directly to simple visibility data. In due 
course, we hope to further develop these simulations to more realistic scenarios and to properly 
compare the technique applied in both image and visibility plane, before eventually applying to 
real data. 
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